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Abstract 


We measured the dipole of the diffuse y-ray background (DGB), identifying a highly significant time- 
independent signal coincidental with that of the Pierre Auger UHECR. The DGB dipole is determined from flux 
maps in narrow energy bands constructed from 13 yr of observations by the Large Area Telescope (LAT) of the 
Fermi satellite. The y-ray maps were clipped iteratively of sources and foregrounds similar to that done for the 
cosmic infrared background. The clipped narrow energy band maps were then assembled into one broad energy 
map out to the given energy starting at E=2.74 GeV, where the LAT beam falls below the sky’s pixel 
resolution. Next we consider cuts in Galactic latitude and longitude to probe residual foreground 
contaminations from the Galactic plane and center. In the broad energy range 2.74 < E< 115.5 GeV, the 
measured dipoles are stable with respect to the various Galactic cuts, consistent with an extragalactic origin. 
The ^-ray sky's dipole/monopole ratio is much greater than that expected from the DGB clustering component 
and the Compton-Getting effect origin with reasonable velocities. At —:(6.5—7)96 it is similar to the Pierre 
Auger UHECRs with Eyyscr = 8 EeV, pointing to a common origin of the two dipoles. However, the DGB 
flux associated with the found DGB dipole reaches parity with that of the UHECR around Eugugcg < 1 EeV, 
perhaps arguing for a non-cascading mechanism if the DGB dipole were to come from the higher-energy 
UHECRs. The signal-to-noise ratio of the DGB dipole is largest in the 5-30 GeV range, possibly suggesting 


the 7-photons at these energies are the ones related to cosmic rays. 


Unified Astronomy Thesaurus concepts: High-energy cosmic radiation (731); Ultra-high-energy cosmic 


radiation (1733) 


1. Introduction 


The energy spectrum of the diffuse 7-ray background (DGB) is 
known accurately from Fermi Large Area Telescope (LAT) meas- 
urements: Fpag = E?dN /dE = 4 x 10-( 7 GeV cm s-! sr-! 
at 0.1S E(GeV) < 100 with 420.35 (Ackermann et al. 
2012b, 2015). The power of its angular anisotropies have 
been measured to be flat (white noise) at Ci Faas £X 
(8-9) x 10™6sr out to 0-7? (£750; Ackermann et al. 
2012a, 2018). The DGB dipole has not yet been measured 
although it would carry important cosmological information arising 
from one of the following three origins: (1) The Compton-Getting 
(Compton & Getting 1935) effect from the Sun's motion generates 
an amplified dipole over that of the cosmic microwave back- 
ground (CMB), assuming the latter is entirely kinematic at 
Vems = 370 kms ^! (Kogut et al. 1993). Since E~?dN/dE is 
Lorentz invariant (Peebles & Wilkinson 1968), the dimensionless 
kinematic DGB dipole would be Dpgp = (4 + p ~ 


0.56%(V / Vemp), an amplification that is comparable (Maoz 1994) 
to that for cosmic rays (CRs; KachelrieB & Serpico 2006) and the 
cosmic infrared background (CIB; Kashlinsky 2005a; Kashlinsky 
& Atrio-Barandela 2022). (2) The clustering (measured to have 
white noise angular spectrum) component dipole from the 
unresolved DGB sources would be E Cr / Fpag œ 0.25% using 
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measurements from Ackermann et al. (2012b, 2018) and also here. 
(3) The Fermi-LAT photons may trace other cosmic structures, 
such as those implied by the UHECR Pierre Auger dipole of 
~(6.5-7)%F pcp (Pierre Auger Collaboration et al. 2017; Aab 
et al. 20202). 

In proper interpretation, the DGB dipole must be described 
not only with amplitude but also with its direction. Here, 
however, fine measurements at high signal-to-noise ratio (S/N) 
are required since the directional uncertainty in the limit of 
S/N>1 is AO c J2(S/N)"! rad (Fixsen & Kash- 
linsky 2011), and only CMB dipole, probed at S/N > 200 
(Kogut et al. 1993; Fixsen et al. 1994), currently allows high- 
precision directional determination. 

We present the first measurement of the source-subtracted 
flux dipole over [3-100] GeV using Fermi 13 yr data. We 
assemble the UltraClean Fermi data set using HEALPix with 
Nsiae = 128 (Górski et al. 2005) from the 13yr all-sky 
observations at narrow energy bands. After removing known 
sources, each E-map is cleaned of the ecliptic plane, the 
Galactic plane, and remaining Galactic and extragalactic 
sources, structured foregrounds, and noise excursions via a 
standard iterative procedure used in CIB work (see, e.g., 
reviews by Kashlinsky 2005a; Kashlinsky et al. 2018). Then, 
possible remaining Galaxy emissions are further removed with 
progressive Galactic cuts in (Iga, baai). The clipped, narrow £E- 
maps were assembled into broad bands of increasing upper 
energy and the residual map dipoles evaluated for each 
situation. We find a robust remaining source-subtracted dipole 
that is highly statistically significant and independent of further 
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Galaxy and source removal, consistent with an extragalactic 
origin. The signal appears stable in time for at least 6 yr. The 
relative amplitude of the dipole is (6—7)96, which is much 
higher than could come from the DGB source clustering 
component and also higher than any kinematic Compton- 
Getting component, unless the local velocity is V 
3000 km s^. However, the relative dipole amplitude coincides 
with that of the UHECR dipole from Pierre Auger Observatory 
(Pierre Auger Collaboration et al. 2017; Aab et al. 20202), 
which we interpret as indicative of a common origin of the 
[3-100] GeV photons here and the UHECR. 


2. Fermi-LAT Data and Processing 


The Fermi-LAT is a pair-conversion telescope covering the 
nominal [0.1—1000] GeV energy range (Atwood et al. 2009). It 
is operated as a sky-survey instrument, having both a wide field 
of view, ~2 sr, and large effective area: ~0.7 m? for ~1 GeV 
photons at near normal incidence. The telescope is comprised 
of a 4 x 4 assembly of modules that each consists of a tracker 
and calorimeter to record information needed to reconstruct 
photon direction and energy. The instrument is enclosed by an 
active anticoincidence system that allows the onboard electro- 
nics to reject charged particle events. 

We selected data taken at the start of scientific operations, 
continuing through the end of mission Cycle 13, specifically, 
all the weekly photon files 2008 August 31—2021 September 1. 
We used the weekly photon files resulting from the latest 
reprocessing, P8R3’ and the front-plus-back converting events. 
To limit contamination from ^-rays scattered from the Earth's 
atmosphere, a selection cut was applied to remove data taken 
when the LAT boresight rocked to >52° with respect to the 
zenith. The selection cuts resulting from this reprocessing 
significantly reduce the occurrence of CR-induced spurious 
(i.e., non-photon) events and importantly, remove residual- 
associated anisotropies (Bruel et al. 2018). It should also be 
noted that these spurious events can introduce anisotropies in 
large- or all-sky analyses. This is due to both noninteracting 
heavy ions and to CR electrons leaking through the ribbons of 
the anticoincidence detector, the latter source being responsible 
for the background anisotropy. We also excluded any time 
interval when the LAT was not in survey mode. To minimize 
potential contamination from CR-induced events in the 
detector, we selected the “UltraClean” event selection cut and 
the corresponding instrument-response function. 

The possible effects on the measured DGB dipole by 
emissions emanating from within the solar system are 
demonstrated below to be negligible. 

Our analysis benefits from including as many celestial 
photons as possible, but bright sources above our latitude cuts 
that are in relatively close proximity would likely introduce 
anisotropies. Thus we exclude emission from the remaining 
brightest 30004- sources from the fourth Fermi Gamma-Ray 
Lat (4FGL) catalog (Abdollahi et al. 2022), then employing the 
clipping procedure as detailed below. 

The 1—100 GeV energy ^-ray sky at low Galactic latitudes is 
dominated by diffuse emissions from interactions of CRs with 
the interstellar medium as well as with the Milky Way radiation 
fields. While these foreground emissions, which have been 
studied extensively (Ackermann et al. 2012b), provide a useful 
diagnostic tool to study the interstellar medium and CR 
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propagation, they are an impediment to large-scale extragalac- 
tic studies. The Galactic latitude bgay and longitude lga cuts 
employed in our analysis partially remove the effects of these 
structured foregrounds on our analysis. We assess further 
effects by making incremental cuts at |bga| > 20°, 30°, 45? and 
30° and at 45? < lca « 330°, 315°, and also by removing the 
ecliptic plane and then recomputing the dipole. 

We apply the clipping procedure well known and commonly 
used in the source-subtracted CIB studies since the COBE/ 
DIRBE, later refined for and further tested in the various 
2MASS, Spitzer, and Euclid work (see reviews in 
Kashlinsky 2005a; Kashlinsky et al. 2018). The method 
isolates iteratively the pixels with photons exceeding a given 
threshold of Npnotcat, bea) 2 (Nphot) F Neuro (N phot)» removing 
here the entire beam at 95% c.l. around the pixels identified in 
the given iteration and proceeds until no more such excursions 
are found at the given Noy. Typically up to ten iterations were 
needed to converge, and the results are insensitive at Nout < 4. 
Of course, the method should be applied to as narrow bands as 
possible, especially for the y-ray sky where adjacent energies 
often trace different sources. Figures 5 and 6 show the 
commonality of the remaining photons for the unclipped and 
clipped maps, demonstrating the efficiency and necessity of 
clipping. After clipping, most remaining pixels starting already 
around (3-5) GeV trace different sources at two adjacent 
bands, and most of the sky becomes progressively darker (no 
photons in each pixel) at narrow energy bands, as expected. 
Next we combine different bands into one broader band map 
and then evaluate the dipoles and the mean remaining flux. 

Figure 1(a) shows the number of photons left in each band 
after the procedures. The horizontal green line in Figure 1(a) 
shows the pixel area (right vertical axis), and the green stars 
mark the 9596 LAT beam area that falls below our pixel 
resolution at E > 3 GeV. Those remaining after clipping y-ray 
sources are unresolved at this resolution. 

We therefore restrict our dipole DGB analysis to 
E> 2/74 GeV. The clipped maps at narrow energy bands are 
then coadded up to E= 115 GeV. Figure 1(b) shows the 
consistency of the histograms of the assembled source-subtracted 
maps for various |bgai|-cuts and 30? < Iga: < 330°. The resultant 
clipped map assembled over 2.74 < E(GeV) < 115 is shown in 
Figure 2; the symbols are explained below. 


3. Dipole Results 


Computations of the multipole expansion of source-sub- 
tracted, band-added maps were done with standard HEALPix 
routines: dipole with REMOVE DIPOLE and power spectrum 
with ANAFAST developed with singular-value decomposition 
(Górski et al. 2005). 

The dipole from the clipped maps is shown with circles in 
Figure 3. The dashed lines of three colors in Figure 1(b) show the 
mean flux in the clipped maps of three values of beut, and the solid 
lines of the same color show the flux dispersion of each map. The 
flux dispersion is or^ 0.8 x 1077 GeV cm? s7! sr |, which 
introduces random dipole uncertainty in each of the three 
components over a map of Ni, pixels of ~op / Nyix. For the 
dipole amplitude, the uncertainty will be v3or~ 5x 
10710 /N../10* GeV em ^s ^ sr", which shows that the 
dipoles found in Figure 3 are statistically significant. A more 
accurate estimate of the statistical error can be done by bootstrap. 
We generated 1000 simulated maps by assigning to each point on 
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Figure 1. (a) Total number of photons in each narrow band is marked with filled circles. Black, blue, and red circles correspond to |bgai| > 20°, 30°, and 45? while 
filled (open) circles corresponds to 0 < Iga < 360° (30? < Ig, < 330°). The LAT beam is shown with green stars for the values displayed on the right vertical axis. 
The horizontal thick green solid line marks the pixel area for N;i4, = 128. The slope of the DGB is marked with black solid line. (b) Flux histograms for 2.74 < 
E(GeV) < 115 and 30? < lca < 330° with black/blue/red corresponding to |bgai| > 20°, 30°, and 45°. Dashed and solid vertical lines of three colors mark the mean 
flux and flux dispersion of each map. 


2.74 € E(GeV) € 115.5 


O.O -—*— 908-07 


Figure 2. The source-subtracted LAT maps coadded from 2.74 to 115 GeV. The displayed sky is masked for |bgai| € 20°. The coordinate grid is marked every 15°. 
The DGB dipole directions are plotted for the 45? < Iga < 315°, |ba| > 45° with red colors and the 30? < Iga < 330°, |bau| > 30° with green colors. The two 
contours in each case correspond to 1c and 20 (68% and 95%) dipole uncertainties. The size of each contour decreases as it should when larger sky area is considered. 
The smaller/narrower contours indicate the directions and uncertainties of the 13 yr assembled map shown here. The larger/broader contours show the same once the 
(A—B)/ 2 dipoles have been subtracted. The yellow arrow indicates the tentative position of the Pierre Auger dipole prior to the deflection by the modeled Galactic 
magnetic field. The Pierre Auger UHECR dipole is marked with the yellow spot and its 68% directional uncertainty is ~15°. The map with the ecliptic plane masked 
out with |8ga] € 5° is shown in Figure 7. 


the unmasked regions of the sky a value extracted at random from Next, we computed the statistical uncertainty in the dipole 
the flux map of 13 yr of data. All data points were given the same direction using the errors c; by generating k= 10* Gaussian 
probability. Following bootstrap standard practice, the selected distributed vectors d* about the dipole components D; measured 
values were not removed. Once a random map j was generated, its from the 13 yr flux map, dé = D; + r*o; with rf Gaussian 
dipole d was computed using the HEALPix remove dipole random numbers with zero mean and unit variance. These 
routine. The mean of each dipole component is zero, and its bootstrap error estimates do not probe the time-variable 
dispersion is a measure of the statistical error. The ratio of the component of the dipole that might be present in the data 
signal to its statistical uncertainty is largest in the energy range from instrument noise, systematics, or genuine y-ray 
~(5-30) GeV. variability. 
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Figure 3. Black, blue, and red colors mark |bg4| 2 20°, 30°, 45? for 
30? € Iga < 330°. The points are slightly shifted at each E for clearer display, 
and the configuration of 45? < Iga, < 315? is not shown to avoid clutter (it is 
shown in Figure 4). The square root of the clustering mean power over £ > 20 
from Figure 3 are shown with squares. The dipoles in the source-subtracted 
13 yr maps are shown with filled circles with the 68% uncertainties. 
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Figure 4. Dimensionless dipole in the same notations as in Figure 3. In 
addition to the numbers in Figure 3, red/blue/black open circles denote the 
configurations with 45° < lca < 315° and 2.74 < E(GeV) < 115 slightly 
shifted in E for clear viewing (this is not shown in the earlier figure to avoid 
overcrowding—all the numbers are consistent with what is shown here). The 
dimensionless dipoles from CRs (UHECR) at 6.5%, Compton-Getting 
amplified CMB dipole, and clustering component are marked with horizontal 
lines. Removing the ecliptic plane, as discussed in the Appendix, further 
reduces the scatter—see Figure 8. The statistical uncertainties are shown at 
95% c.l. 


For estimating the time variations and noise in the maps, we 
assembled two time-separated subsets of data: A for years 1—6 
and B for years 7-12 with the odd year 13 left out temporarily 
for this task. Each of the two 6 yr time intervals are longer than 
the typical variability timescales of y-ray sources. The masks 
from the clipping of sources and Galaxy for the full 13 yr map in 
(baa, laa) are applied to each subset, and the narrow E-bands are 
then coadded to produce the 13 wide band maps going from 
Eo = 2.74 GeV to the given E up to E = 115 GeV. Then, similar 
to the by-now standard CMB (Smoot et al. 1992; Bennett et al. 
1996) and CIB (Arendt et al. 2010; Kashlinsky et al. 2012) noise 
power processing, the noise properties are computed from the 
(A — B)/2 maps, masked in the same way as the combined 13 yr 
data. The resultant noise map is shown in Figure 9, and the 
dipole amplitudes of the 13 yr and (A — B)/2 and dipole powers 
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are compared in Figure 11. Figure 10, left panel, displays the 
histograms in the time-differenced maps with the ecliptic plane 
masked out, showing the standard deviations using 12 yr data of 
Cr-By2 ~ 1.2 x 10" GeVcm ?s !sr !. The right panel 
plots the correlation coefficient between the (A — B) and the 
final 13 yr maps and demonstrates negligible-to-null correlations, 
implying that the two maps are independent with the dipole 
power from the (A — B)/2 contributing in quadrature to the 
signal shown in Figure 3. 

Regarding the possibility of systematics, our analysis is not 
dependent on a precision spectral determination or on an 
absolute flux, both of which are susceptible to uncertainties in 
the instrument response. Furthermore, we are not subtracting a 
Galactic foreground model, which may introduce potential 
systematics. Finally, our error analysis based on the A — B time 
differences argues against the presence of time dependencies 
systematics, which could arise from systematic effects in the 
exposure maps. 

We evaluated the power spectrum of the dipole- and 
monopole-removed source-subtracted 13 yr maps and find it 
is approximately white out to the computed /— 500, in 
agreement with the earlier analyses (Ackermann et al. 
2012b, 2018), except about 20%-25% lower because of our 
(more efficient) source subtraction. The clustering dipole 
contribution, from averaging the (white) power at 
20 < £ < 500 to reduce the effects of the cosmic variance at 
low £, is shown with squares in Figure 3. For the clustering 
independent of the uncovered dipole, the latter would add in 
quadrature contributing negligibly in the final dipole balance 
compared to the found dipole. 

The angular uncertainty on the measured y-ray dipole is shown 
in Figure 2 for two configurations |ba4| = 30°, 30? < Iga < 330° 
(green contours), and |bg4| — 45?, 45° < Iga < 315? (red con- 
tours). To estimate the directional uncertainty, we likewise 
generated k— 10^ random dipoles d' with the uncertainties 
derived from bootstrap. The angular separations between the 
measured and simulated dipoles sorted in increasing order and the 
lo, 2c confidence intervals were defined as the regions that 
encloses 68% and 95% of all angular separations. The 
uncertainties are marked by the two sets (red and green) of small 
contours. This dipole determination includes diffuse and time- 
variable contributions. To estimate the direction and its 
uncertainty of the diffuse y-ray component alone, we subtracted 
the dipole vector D(4—2)/2 from the 13 yr dipole vector D13 yr and 
added the errors in quadrature. We repeated the procedure as 
before, and the results are shown by the wider red and green 
contours. Although the time-varying y-ray component on the 
direction has some effect, in the Appendix we show that its effect 
on the dipole amplitude is negligible. The Pierre Auger dipole is 
shown as a yellow point, and the arrow indicates its possible 
location prior to the deflection produced by a model of the 
Galactic magnetic field. 


4. Discussion 


Figure 4 summarizes our findings in terms of the 
dimensionless dipole amplitudes. The identified dipole has 
dimensionless amplitude ~(6.5—7)% out to E ~ 100 GeV and 
can be either of Galactic origin or extragalactic from one of the 
three origins mentioned in Section |. It is unlikely to come 
from the Galaxy as its properties show no dependence on the 
Galactic cuts, |bgail, lGa outside of the Galactic center for 
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E = 3 GeV. Foregrounds are removed efficiently here, and the 
DGB is reproduced for the remaining sky as Figure 2 shows. 

The dipole directions are consistent for all Galaxy cuts and E 
configurations with the 6896 c.l. errors of ~10°-20°. The rms 
directional uncertainty is close to AO ~ 4/2 (S/N)^! radian 
(Fixsen & Kashlinsky 2011). For data in the region 
30? < Iga < 330° and cuts in |bgai| = 20°, 30°, 45°, the dipoles 
are located in the range (lca, baa) = (—30° — ~ 0); for the 
data in the region 45° < / < 315° and for the same cuts in |bga1| 
the dipoles are at Galactic (J, b) = (270° — 300°, — 40? — 0). 
In Figure 2, we show the direction and its 68% and 95% 
angular dispersion for two configurations with 30? < | < 330°: 
cuts at |bgai| = 30° and |bgai| = 45°. We note that the direction 
comparison is further hampered by the deflection of CRs by the 
Galactic magnetic field, shown by the (model-dependent) 
location of the UHECR dipole prior to this effect. 

We tested for solar system effects, primarily from the Sun 
and Moon, by removing the ecliptic plane and find no 
noticeable changes in the dipole as discussed in the Appendix. 
The distribution of solar system emissions are expected to be 
sharply distributed about the ecliptic (Johannesson & 
Orlando 2013). We considered the possibility of y-ray emission 
emanating from the solar system by masking the 10° band of 
the ecliptic plane and find the same signal within the errors, 
indicating negligible contribution from there demonstrating at 
most a negligible dipole contribution to the measured signal. 

The expected dipole component from clustering, detustering, 1S 
shown in Figure 3 to be on average ~0.2% from using the power 
spectra computed here, which is consistent with Ackermann 
et al. (2012b, 2018); if the found dipole is independent of 
clustering, deustering is added in quadrature, making this 
contribution totally negligible. Cosmic variance affects the 
measured clustering dipole per probability distribution 


Sustering @Sustering : 
P(detustering) = 4 oe ew[ J (Kaiser 1983; Abbott & 


26 G 26 


Wise 1984); hence, at 9596, 99% c.l., detustering < 4, 5746. 
which is still well below the found dipole. This makes clustering 
an unlikely origin for the measured dipole. The fact that the 
measured diffuse emission power spectrum agrees well with the 
other measurements and shows no drastic dependence on the 
Galaxy cuts further indicates that the diffuse emission, and likely 
its dipole, do not have critical Galactic components. 

If the dipole originates from the Compton-Getting effect due 
to our kinematic motion, the implied velocity would be 
2:3000 km s^, an order of magnitude larger than the CMB 
dipole equivalent velocity marked with the dashed line in 
Figure 4 and significantly larger than the effective velocities 
implied in other dipole measurements (Kashlinsky et al. 
2008, 2010, 2011; Singal 2011; Atrio-Barandela et al. 2015; 
Secrest et al. 2022). Thus, the Compton-Getting amplified 
kinematic origin of the measured DGB dipole also appears to 
us implausible. 

At face value, there appears a commonality with the 
UHECR dipole from the Pierre Auger Observatory (Pierre 
Auger Collaboration et al. 2017). The dimensionless ampl- 
itude of the UHECR dipole at E = 3 EeV is shown in Figure 4 
with the solid line at 6.5%. It coincides well with the found 
DGB dipole, perhaps suggesting a common origin. The source 
of the UHECRs in general is still unknown. Active galactic 
nuclei (AGNs) have long been considered as potential sites 
for production of high-energy CRs, including possibly 
UHECRs (Fang & Murase 2018; Murase & Stecker 2022). 
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Also, there has been similar speculation regarding y-ray bursts 
(GRBs), and they cannot currently be ruled out as sources of 
prolific UHECR emission (Murase et al. 2022). However, 
both of these populations are well known to be highly 
isotropic. So, some as-yet-unidentified subpopulation would 
be required to produce the dipole anisotropy documented by 
the Pierre Auger experiment if AGNs or GRBs are its source. 
Very generally, a common origin of UHECR and the 
[3-100] GeV DGB may be due to pionic photons (Halzen & 
Hooper 2002) from photomesonic production (Stecker & 
Salamon 1999) arising in several mechanisms. Possible 
venues may come from UHECRs above the Greisen- 
Zatsepin-Kuz min (GZK) knee around ~10°° EeV (Grei- 
sen 1966; Zatsepin & Kuz’min 1966), with pions arising from 
the A resonance: p+ cmp — At — q^ +p or more spec- 
ulatively from proton decay (Sakharov 1967) p — e* + 7°, 
pop + 7° (Tanabashi et al. 2018). Subsequently, no 2y, 
and on average, the pionic photons would carry ~1/2 of the 
parent pion’s energy (Mészáros 2014; Halzen & Kheirand- 
ish 2019; Globus & Blandford 2023). GZK horizon of <100 
Mpc would then apply to the sources in the latter case 
(Stecker 1968; Ding et al. 2021) whose precise value depends 
on the CR composition (Ahlers & Salvado 2011). The pions 
trigger electromagnetic cascades, transferring energy from the 
CRs to the GeV-TeV y-ray photons (Fornasa & Sánchez- 
Conde 2015), but their energy flux would have to be up to 
50% of the CR flux (Kalashev et al. 2009). However, the 
dipole here corresponds to the [3-100] flux of 
~(1-1.5) x 107 GeV cm 7s 'sr-', which is reached for 
UHECRs around Eyyecr S 1 EeV (Mészáros 2014; Fang & 
Murase 2018; Aab et al. 2020b), so it is not clear whether a 
cascading mechanism can produce the y-ray signal out of the 
UHECR protons. More likely this points to a common origin 
between the CRs and the high-F photons observed here by 
Fermi-LAT. See possible specific models in, e.g., Waxman 
(1995); Fang & Murase (2018); Ding et al. (2021). In this 
case, the distance to the sources may not exceed the horizon 
for high-E y-rays (Nikishov 1962; Stecker 1968; Fazio & 
Stecker 1970) from the CIB and optical background 
(Kashlinsky 2005b; Helgason & Kashlinsky 2012; Fermi- 
LAT Collaboration et al. 2018). These sources lie at smaller 
distances than the cluster sample used in the dark flow dipolar 
probe (Kashlinsky et al. 2008, 2010) or the source sample of 
the WISE dipole probe (Secrest et al. 2022) and cannot affect 
those analyses. 

However, as a result of magnetic deflection, using UHECR 
anisotropy information to constrain the spatial distribution or 
the nature of the CR sources is not plausible given the 
limitations inherent in galactic and extragalactic magnetic field 
models (Allard et al. 2022). Indeed, we cannot unambiguously 
claim an association between the UHECR and ^-ray dipoles 
although the coincidence in amplitude and nominal orientation 
is intriguing. We would note, though, that in the case of a 
common source, the dipole nature of the UHCER flux is likely 
to be preserved despite magnetic deflection of individual 
particles. 

In addition to the issue of magnetic deflection it must be 
considered that the sources of ^-rays and UHECRs may not be 
strictly in common. For example, it is apparently the case that 
some subset of ^-ray-emitting blazar AGNs are prolific sources 
of high-energy neutrinos, while some are not; or perhaps the 
production sites within the latter are opaque to y-ray emission 
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(see, e.g., Plavin et al. 2023). Additionally, there are apparently 
non-blazar sources of high-energy neutrinos, e.g., NGC 1068 
and perhaps the Milky Way itself. The Ice-Cube neutrinos may 
be unrelated to the Auger UHECRs, but by analogy, the 
underlying UHCER sources may comprise a_ similarly 
nonhomogeneous population. 
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Appendix 


A.1. Map Clipping 


Figure 5 shows the commonality of photons in unclipped 
maps at the given E between each two adjacent bands starting 
at 1.15 GeV for the energy range at which the DGB dipole will 
be assessed. The red/magenta color corresponds to the most 
photons, with blue corresponding to the fewest. The figure 
shows that as one reaches higher E, there are fewer and fewer 
common pixels even among the two most adjacent bands. 
Figure 6 shows the same after clipping. 


A.2. Testing for Null Systematics from the Plane of the Ecliptic 


In principle, the Sun produces diffuse 4-rays from the 
inverse Compton scattering on its radiation field, and the Moon 
produces some diffuse y-ray component from CR interactions. 
Once averaged over the year, the Sun and the Moon may 
contribute some anisotropic diffuse components contributions, 
which we test here. 

From, e.g., Ackermann et al. (2015), it is evident that the 
solar system flux can be at a level of up to ~5% of the isotropic 
background level. However, this emission, which is due to CR 
interactions with the solar disk and the Moon as well as the 
inverse-Compton radiation involving the solar radiation field, is 
concentrated along the ecliptic plane; see, e.g., Figure 2, 
Johannesson & Orlando (2013). We also generated emission 
templates as described therein and examined the latitude 
profiles. Examination of the resulting latitude profiles suggests 
that a 10° ecliptic latitude cut should lead to eliminating >90% 
of the solar system emissions. In any case, the ecliptic latitude 
cut did not lead to a discernible difference in our dipole 
determination as noted in Section 4. 
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Figure 7 shows the same map as used in Sections 2 and 3, 
Figure 3, demonstrating no obvious contributions from 
emissions within the ecliptic plane. 

Figure 8 shows the dipole/monopole ratio, as in Figure 4, 
when emissions from the ecliptic plane are not masked out. The 
two sets of numbers are essentially identical within the errors, 
demonstrating no substantial contributions to the identified 
DGB dipole from the solar system and Moon ^-ray emissions. 


A.3. Testing for Time Variability with Time- differenced Maps 


For estimating the time variations and noise in the maps, we 
assembled two time-separated subsets of data: subset A for 
years 1—6 and subset B for years 7-12 with the odd year 13 left 
out temporarily in this task. The two 6 yr time intervals are each 
longer than the variability timescales of y-ray sources. The 
masks from the clipping of sources and Galaxy for the full 
13 yr map in (bca, lca) are applied to each subset, and the 
narrow E-bands are then coadded to produce the 13 band maps 
going from Eo — 2.74 GeV to the given E up to the final 
E — 115 GeV. 

Figure 9 shows the resultant time-differenced map, (A — B)/2, 
with the Galactic plane cut out at |bo4] = 20? 

Figure 10, left panel, displays the histograms in the time- 
differenced maps with the ecliptic plane masked out. The 
source-clipping mask from the full 13 yr data sets has been 
applied, and the sky is kept at 30? < Iga < 330? and bg 2 20°, 
30°, 45° for the black/blue/red colors. The standard deviations 
are shown with vertical solid lines corresponding to 
OFA—B)/2~1.2 x 1077 GeV cm ? s! sr '. 

Figure 10, right panel, plots the correlation coefficient, 
Ra gy x 13 yr, between the (A—B) and the final 13 yr maps. The 
figure demonstrates negligible-to-null correlations, implying 
that the two maps are independent with the dipole power from 
the (A—B)/2 contributing in quadrature to the signal shown in 
Figure 3. The precise amplitude, and sign, of the negligible 
cross correlation appears sensitive to the (small) addition of ^- 
ray photons at higher E. 

Figure 11 shows the ratio of the dipole powers, G = (d), in 
the (A—B)/2 to that in the 13 yr maps. With the uncorrelated 
time-varying contribution to the 13yr dipole added in 
quadrature, this implies that the dipole power contributed, C4, 
from the time-differenced maps is well below < 10%. Whereas 
the 13 yr maps appear to have a stable dipole for the various 
considered configurations, the (A — B)/2 dipole power varies 
drastically, by more than 2 orders of magnitude, further 
suggesting the absence of their substantive contributions to the 
uncovered dipole in Figure 3. 

Thus, the time variations between the two 6 yr subsets 
cannot make appreciable contribution to the dipole power 
found here. 
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Figure 5. Number of photons per given pixel common in the two adjacent maps is shown with density plots for unclipped maps. The pixels where there are zero 
photons in one of the two bands are shown at n, = 0.01. The shown colors span 0-5500 common pixels logarithmically. 
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Figure 6. Same as Figure 5 but for clipped maps. 
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Figure 7. Map with the ecliptic plane masked as in Figure 4 without masking the ecliptic plane. 
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Figure 8. Dimensionless dipole from the maps in Figure 7 in the same notations as in Figure 4 without masking the ecliptic plane. The error bars are barely 
distinguishable from those in Figure 4 and are not shown to avoid clutter. 
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2.74 < E(GeV) < 115.5, (A-B)/2 
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Figure 9. The (A — B)/2 map from 2.74 to 115 GeV is shown with the overall 13 yr clipping mask applied and Galaxy removed for |beai| < 20°. 
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Figure 10. (Left) Histograms of fluxes after the 13 yr clipping mask has been applied to the maps shown in Figure 9. The standard deviations are shown with vertical 
solid lines corresponding to oF (4-8)/2~1.2 x 107 GeV cm ^s ' sr. The masking from clipping has been applied, and the sky has been kept at 30° < Iga, < 330° 
and |bgai| > 20°, 30°, 45° for the black, blue, and red colors. (Right) The correlation coefficient, R, between the (A — B) and the final 13 yr maps. Black, red, and blue 


symbols correspond to |b.| > 20°, 30°, 45°. The filled and open symbols correspond to 30° < /54 < 330° and 45° € Iga € 315°, respectively. 
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Figure 11. The ratio of the dipole power, Cj = (dẹ), in the (A — B)/2 to 13 yr 
maps. Open circles correspond to 45° < lca x 315? and filled to 
30° < Iga < 330°. Black, blue, and red denote |bgqi| > 20°, 30°, 45°, 
respectively. 
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